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SECTION  I 


*  INTRODUCTION 

1.  Background 

The  determination  of  the  structural  response  induced  by  thermal  effects 
is  an  important  factor  in  many  aerospace  structural  designs.  Extreme  aero¬ 
dynamic  heating  on  advanced  aerospace  vehicles  may  produce  severe  thermal 
stresses  that  can  reduce  operational  performance  or  even  damage  structures. 
The  performance  of  laser  devices  can  be  degraded  by  thermal  distortions  of 
mirror  surfaces.  The  thermal  environment  in  space  may  cause  orbiting 
structures  to  distort  beyond  operational  tolerances.  To  predict  the 
structural  response  accurately,  effective  nimterical  techniques  capable  of 
both  thermal  and  structural  analyses  are  required.  One  technique,  the 
finite  element  method,  has  been  found  to  be  particularly  well-suited  for 
such  analyses  due  to  its  capability  to  model  complex  geometry  and  to  perform 
both  thermal  and  structural  analyses. 

In  predicting  the  thermal-structural  response,  basic  thermal  elements 
with  assumed  linear  temperature  distribution  are  frequently  employed.  Nodal 
temperatures  obtained  from  the  thermal  finite  element  analysis  are  trans¬ 
ferred  to  the  structural  finite  element  analysis  for  computations  of  dis¬ 
placements  and  stresses  where  elements  with  linear  displacement  distribu¬ 
tions  are  used.  This  procedure,  denoted  as  the  conventional  approach,  is 
shown  schematically  in  Figure  1(a).  With  the  use  of  the  elements  with 
linear  distributions  in  both  thermal  and  structural  analysis,  a  large  nunber 
of  el ements  are  normally  needed  to  produce  accurate  thermal -structural  solu¬ 
tions.  Often,  basic  differences  in  the  thermal  and  structural  problems  man¬ 
date  different  analysis  models,  and  the  data  transfer  between  the  analyses 
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can  become  complicated  because  of  the  need  to  interpolate  temperatures  at 
the  structural  nodes.  With  increasing  structural  design  complexity  and  the 
need  for  highly  accurate  analysis,  improvements  in  finite  element  methods 
are  needed  to  increase  the  accuracy  and  efficiency  of  coupled  thermal-struc¬ 
tural  analysis. 

2.  Objectives 

To  improve  the  accuracy  and  efficiency  of  the  finite  element  method, 
development  of  an  approach  called  integrated  thermal-structural  analysis  was 
initiated  in  References  1-5.  The  goals  of  the  integrated  approach  are  to: 
(1)  provide  thermal  elements  which  predict  detailed  temperature  distribu¬ 
tions  accurately,  (2)  provide  structural  elements  with  improved  displacement 
and  stress  distributions  which  are  fully  compatible  with  the  thermal  ele¬ 
ments,  and  (3)  integrate  the  thermal  loads  with  the  structural  analysis  to 
further  improve  the  accuracy  of  displacements  and  stresses.  These  concepts 
are  shown  schematically  in  Figures  1(b)  and  1(c). 

The  goals  of  the  integrated  approach  require  developing  new  thermal  and 
structural  finite  elements  that  can  provide  higher  accuracy  and  efficiency 
than  conventional  finite  elements.  The  task  of  developing  new  thermal  ele¬ 
ments  with  a  common  structural  element  discretization  was  begun  in  Refer¬ 
ences  1-5.  Displacements  and  stresses  based  on  the  new  thermal  elements 
were  improved  because  more  accurate  thermal  loads  were  provided  to  the 
structural  finite  element  analysis. 

Another  task  of  the  integrated  approach  is  to  develop  new  structural 
elements  capable  of  providing  improved  displacement  and  stress  distribu¬ 
tions.  By  integrating  these  new  structural  elements  with  the  new  thermal 
elements  developed  previously,  better  thermal-structural  solution  accuracy 
can  be  obtained. 
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Figure  1.  Method  of  improving  thermal -structural  solution  accuracy. 


3.  Scope 

To  further  develop  integrated  thermal-structural  analysis,  two  proce¬ 
dures  for  developing  improved  thermal-structural  analysis  are  presented.  A 
brief  review  of  the  conventional  thermal-structural  finite  element  formu¬ 
lations  is  first  given  in  Section  II.  This  section  describes  element  chat - 
acteristics  including  some  deficiencies  which  motivated  the  development  of 
improved  structural  elements.  The  new  approaches  for  improving  the  accuracy 
and  efficiency  of  thermal -structural  finite  elements  are  presented  in  Sec¬ 
tions  III  and  IV.  An  approach  using  a  hierarchical  finite  element  concept 
is  introduced  in  Section  III,  and  an  approach  using  nodeless  parameter 
structural  finite  elements  is  introduced  in  Section  IV.  In  Section  V,  some 
benefits  of  utilizing  the  new  integrated  thermal-structural  approaches  are 
demonstrated  by  numerical  examples.  Section  VI  discusses  the  results  and 
highlights  areas  for  future  research. 
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SECTION  II 


CONVENTIONAL  THERMAL-STRUCTURAL  FINITE  ELEMENTS 


The  objectives  of  this  section  are  to  briefly  describe  conventional 
thermal-structural  finite  element  formulations  and  demonstrate  the  behavior 
of  thermal  stresses  obtained.  These  basic  finite  element  formulations  will 
be  referred  to  in  later  sections.  For  simplicity  in  understanding  the 
thermal  stress  behavior,  a  rectangular  element  shape  in  two  dimensions  is 
used.  Details  of  finite  element  formulations  for  general  quadril ateral 
element  shapes  appear  in  Reference  5. 

1.  Thermal  Element  Formulation 
1.1  Element  Interpolation  Functions 

The  element  temperature  distribution  for  a  conventional  bilinear  four- 
node  thermal  element  in  local  Cartesian  coordinates  is  expressed  in  the 
form, 


(1) 


where  N^  and  T^,  i  =  1,4  are  the  element  interpolation  functions  and 
the  time  dependent  nodal  temperatures,  respectively.  For  a  rectangular 
element,  the  element  interpolation  functions  are  defined  by, 


N,  -  (1  -  -Ml  -  i) 


n2  -  i  (l  - 1) 


a  b 


a  b 


N4  =  (1  -  — )  — 


a  b 


(2) 


a  b 
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where  a  and  b  are  the  element  dimensions  in  the  x  and  y  directions. 
With  these  element  interpolation  functions,  an  element  temperature  distri¬ 
bution  varies  bilinearly  as  shown  in  Figure  2. 

1 .2  Element  Matrices  and  Heat  Load  Vectors 

Finite  elanent  (F.E.)  formulations  for  nonlinear,  transient  thermal 
problems  can  be  derived  from  the  governing  heat  conduction  equation  with 
radiation  boundary  conditions  by  the  method  of  weighted  residuals  (Reference 
6).  In  general,  the  element  temperature  T(x,y,t)  and  temperature 
gradients  are  expressed  in  the  form 


(3a) 


T  ■  [NT]{T(t)}e 


(3b) 


where  ( T( t ) } e  denotes  a  vector  of  element  nodal  temperatures  as  a  function 
of  time.  For  simplicity,  conduction  with  only  specified  surface  heating  and 
radiation  heat  transfer  will  be  considered.  Finite  element  thermal  analy¬ 
ses  for  other  heat  loads  such  as  internal  heat  generation  and  surface  con¬ 
vection  are  presented  in  References  1-2.  For  transient  thermal  analysis  the 
equations  for  a  typical  element  are 


[C]e{T}e  +  [Kc]eme  ♦  [Kr]e{T}e 


=  I Qql  e  +  {Qr}e 


(4) 


$ 


Figure  2.  Two-dimensional  conventional  bilinear  thermal 
element  interpolation  function. 
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where  the  element  matrices  are  expressed  in  terms  of  integrals  over  an  ele¬ 
ment  volune  V  and  surface  S  .  The  element  equations  are 


[C]e  ■  ;v  0C[NT]T[NT]dV 

0 

(5a) 

[Kc]e  =  V  [BT]V][BT]dV 

(5b) 

L{T.)  -  /s  oer*[NT]TdS 

(5c) 

e 

{Qq]e  -  / Se  q[NT]TdS 

(5d) 

{ Qr 1 e  •  /s  aqr[NT]TdS. 

(5e) 

All  thermal  parameters  may  be  temperature  dependent  in  general,  but  are 
assuned  constant  herein. 

2.  Structural  Element  Formulation 
2.1  Element  Interpolation  Functions 

The  two-dimensional  conventional  bilinear  four-node  structural  element 
has  two  in-plane  di spl acements  u  and  v  which  may  vary  with  the  element 
local  coordinates  x,y  and  time  t.  Element  displacement  distributions  are 
assumed  in  the  form. 
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» 


{6}  - 


uCx.y.t)1 

”  Ni  0  N2  0  N3  0  N4  0 

r  = 

.vU.y.t) 

_0  Ni  0  N2  0  N3  0  Ni+_ 

ul 

Vi 

U2 

v2 

u3 

V3 

U4 

V** 


*  = 


[Ns](«}  (6) 


where  N-,  i  =  1,4  are  the  element  displacement  interpolation  functions 
which  have  the  same  form  as  for  the  conventional  finite  element  temperature 
interpolation  functions  shown  in  Equation  2.  Typical  element  displacement 
distributions  are  shown  in  Figure  3. 

2.2  Element  Stiffness  Matrix  and  Thermal  Force  Vector 

With  the  assumed  element  displacement  distributions  shown  in  Equation 
6,  the  element  strain-displacement  relations  can  be  written  as, 


W  = 


au_ 

Cp 

X 

3x 

9v 

£ 

i  “ 

_ — 

y 

3y 

Y 

3  u  (  9v 

,  xy 

9y  3x 

-  [BS]{«1 


(7) 


where  [B<~]  is  the  strain-displacement  interpolation  matrix.  For  quasi¬ 
static  analysis,  the  principle  of  minimum  potential  energy  is  applied  to 
derive  the  element  equations.  Typical  element  equations  may  be  written  in 
the  form 


[Ks]{S}  =  {Ft} 


(8) 
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Figure  3.  Two-dimensional  conventional  structural  element 
interpolation  functions. 
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where  [1C]  is  the  element  stiffness  matrix,  and  { F-j. }  is  the  equivalent 
nodal  thermal  load  vector.  These  matrices  are  expressed  in  the  form  of 
integrals  over  the  element  volume  V  as 


[Kj]  *  Jv  [Bs]T[D][8s]dV 


(9a) 


(Ft|  =  /v  [Bs]T[D]M(T(x,y,t)-Tref)dV  (90) 


where  [D]  is  the  elasticity  matrix  defined  by  (plane  stress), 


[D] 


E 

l-v2 


1 

v 

0 

. 


v  0 
1  0 

0  ii 

2 


(10) 


where  v  is  Poisson's  ratio.  The  vector  {a}  contains  the  thermal  expan¬ 
sion  coefficients  given  by  (plane  stress) 


W 


(ID 


T(x,y,t)  is  the  element  temperature  distribution  defined  in  Equation  1,  and 

T  ^  is  the  reference  temperature  for  zero  stress.  The  elasticity  matrix 
ref 

[D]  and  the  vector  of  thermal  expansion  coefficients  shown  in  the  above 
equations  can  be  used  for  plane  strain  by  substituting  E/(l-\>2)  for  E, 
v  { 1  -v )  for  v  ,  and  { 1  +v  )a  for  a . 
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2.3  Element  Stresses 


After  the  element  matrices  shown  in  Equation  8  are  assembled  and  the 
element  nodal  displacements  (a}  are  computed,  the  element  stresses  can  be 
obtained  using  the  thermo-elastic  stress-strain  relations, 

■  [D]{[Bs]{S)  -  {“} (T(x,y,t)  -  Tref)  |  .  (12) 

3.  Behavior  of  Element  Thermal  Stresses 

The  behavior  of  the  conventional  element  thermal  stresses  can  be  ex¬ 
plained  using  the  thermo-elastic  stress-strain  relations  shown  in  Equation 
12.  These  element  stresses  may  be  simply  written  in  the  form  of  the  differ¬ 
ence  between  element  total  strains  and  element  thermal  strains  as 


fa 

e  -  a  (T-T  ,) 

x  v  ref 

X 

0 

y 

=  [D]  . 

e  -  a  {T-T  .) 

y  v  ref' 

i 

Y  -  0 

xy 

) 

(13) 


For  the  bilinear  element  displacement  distributions  shown  in  Equation  6  the 
total  strain  ex,  the  derivative  of  the  u-displ acement  with  respect  to  x, 
is  constant  in  the  x-direction.  But  the  thermal  strain  a(T-Tref.)  which 
varies  directly  with  the  element  temperature  distribution  is  linear  in  the 
x-direction.  Because  the  total  strain  and  the  thermal  strain  have  different 
polynomial  orders,  unrealistic  element  thermal  stress  distributions  can 
resul t. 
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To  examine  the  conventional  finite  element  thermal  stress  behavior,  a 
fundamental  two-dimensional  thermal  stress  problem  is  analyzed.  Figure  4(a) 
shows  a  plate  free  to  expand  with  a  prescribed  linear  temperature  distri¬ 
bution.  The  plate  is  modeled  by:  (1)  one  conventional  finite  element,  and 
(2)  four  equal  conventional  finite  elements.  Since  the  plate  is  subjected 
to  a  linear  temperature  distribution  and  is  free  to  expand,  all  stress  com¬ 
ponents  should  be  zero.  Figure  4(b)  shows  the  actual  distributions  of  the 
stress  component  obtained  from  these  two  finite  element  models.  The 
stress  distributions  vary  linearly  within  the  elements  and  are  maximum  at 
the  element  boundaries,  although  the  stresses  at  the  element  centroids  are 
zero.  By  increasing  the  number  of  elements  used,  the  element  stress  distri¬ 
butions  are  improved,  but  the  stress  discontinuities  across  the  element 
boundaries  still  appear.  Conventional  finite  elements  provide  accurate 
stresses  at  the  element  centroids  in  this  simple  example.  However,  the 
element  centroidal  stresses  are  not  always  conservative  estimates  of 
stresses  in  thermal  stress  problems  as  will  be  shown  by  examples  in  Section 
V. 

It  should  be  emphasized  that  the  deficiency  of  the  conventional  finite 
element  in  predicting  unrealistic  element  thermal  stress  distributions  is 
because  the  element  total  strain  contains  a  mechanical  strain  term  one  order 
lower  than  the  thermal  strain  term.  Thus,  to  improve  the  element  stress 
distributions,  the  element  displacement  distributions  should  vary 
quadratical ly  so  that  the  element  total  strain  is  linear  as  the  thermal 
strain. 
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500°F 


(alFree  expansion  plate  with  linear  temperature  distribution. 


(b)  Comparative  thermal  stress  solution. 


Figure  4.  Deficiency  of  conventional  finite  element  in  predicting 
thermal  stress  distribution. 


14 


SECTION  III 


HIERARCHICAL  THERMAL-STRUCTURAL  FINITE  ELEMENTS 
1.  Integrating  Thermal  and  Structural  Analyses 

To  more  fully  develop  the  potential  of  the  finite  element  method  for 
thermal-stress  analysis,  the  concept  of  integrated  thermal-structur al  analy¬ 
sis  was  proposed  in  Reference  1.  The  approach  focuses  on  aerospace  appli¬ 
cations  where  often  thermal  and  structural  models  differ  because  of  differ¬ 
ent  analysis  requirements.  The  objectives  of  the  approach  are  to  provide 
more  efficient  coupl  ing  between  the  thermal  and  structural  analyses  and  to 
improve  the  accuracy  of  each  analysis  particularly  the  thermal-stress  analy¬ 
sis.  The  basic  philosophy  of  the  approach  is  that  optimum  accuracy  and 
efficiency  will  be  achieved  by  working  from  a  common  geometric  model  and 
performing  the  analyses  with  the  finite  element  method.  The  technology  for 
creating  geometric  models  within  computer  aided  design  is  relatively  mature 
and  is  a  natural  approach  for  initiating  many  interdisciplinary  analyses 
including  thermal  and  structural  analyses.  Modern  computer  software  such  as 
the  popular  PATRAN-G  program  are  based  on  the  concept  of  first  developing  a 
geometrical  model  and  then  creating  the  finite  element  discretization.  The 
finite  element  method  is  the  logical  choice  as  the  analysis  method  because 
it  is  inherently  based  on  geometry. 

The  integrated  thermal-structural  approach  advocated  herein  is  based 
upon  using  the  sane  geometric  model  for  both  analyses  with  a  common  finite 
element  nodal  discretization.  Ihe  thermal  and  structural  models  are  per¬ 
mitted  to  differ  because  heat  transfer  and  structural  response  may  depend  on 
different  characteristics  of  the  system.  The  heat  transfer  problem  may 
mandate  modeling  features  such  as  insulation  or  cooling  passages.  The 


15 


structural  response,  however,  will  not  depend  on  such  non-load  bearing 
features  of  the  thermal  model.  Only  in  the  case  of  conduction  heat  transfer 
with  simple  boundary  conditions  can  the  thermal  and  structural  models  act¬ 
ually  be  identical.  Thus  an  integrated  analysis  approach  must  permit  the 
use  of  different  finite  element  models  for  the  thermal  and  structural  analy¬ 
ses. 

A  second  consideration  in  developing  integrated  thermal-structural 
analysis  is  to  recognize  that  temperature  and  stress  distributions  may  be 
significantly  different.  Regions  of  high  temperature  gradients  do  not 
always  correspond  to  regions  of  high  stress  gradients.  Temperatures  may 
vary  (or  not  vary)  significantly  over  a  structure  due  to  radiation  boundary 
conditions  or  convective  cooling.  Stresses  may  vary  significantly  in  a  dif¬ 
ferent  region  of  the  structure  because  of  boundary  constraints  and  geometri¬ 
cal  effects  that  cause  stress  concentrations.  Moreover,  the  locations  of 
the  regions  of  high  thermal  or  stress  gradients  are  not  known  apriori.  Thus 
an  integrated  thermal-structural  analysis  must  have  the  capability  for 
independently  refining  the  thermal  and  structural  analysis  to  capture  these 
gradients.  The  approach  advocated  for  refinements  is  based  on  the  use  of 
finite  element  nodeless  variables  and  hierarchical  interpolation  functions. 
Further  details  of  this  approach  are  given  in  Section  III. 

A  third  important  consideration  in  integrated  thermal-structural  analy¬ 
sis  is  effectively  coupling  the  analyses.  The  details  of  the  temperature 
variations  obtained  in  the  thermal  analysis  should  be  employed  consistently 
in  the  structural  analysis  to  obtain  the  true  distributions  of  displacement 
and  stress.  The  finite  element  method  permits  consistent  coupling  through 
equivalent  nodal  thermal  force  integrals.  The  integrated  thermal-structural 
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analysis  approach  consistently  employs  these  integrals  with  the  temperature 
representations  from  the  thermal  analysis. 

2.  Improving  Finite  Element  Approximations 

Three  basic  approaches  for  improving  the  quality  of  finite  element  ap¬ 
proximations  are  to: 

(1)  introduce  additional  elements  of  a  smaller  size  in  an  area  of 
significant  dependent  variable  variation.  This  approach  is  called 
the  h  method  where  h  denotes  the  element  size,  and  it  is  the 
standard  approach  used  for  the  solution  of  practical  problems  via 
standard  production  codes. 

(2)  use  the  same  size  and  definition  of  elements  but  increase  the 
order  of  the  interpolating  polynomials.  This  approach  is  called 
the  p  method  where  p  denotes  the  order  of  the  element  poly¬ 
nomials.  With  this  approach  additional  degrees  of  freedom  are 
introduced  via  nodeless  variables.  The  interpolating  polynomials 
progress  in  order  from  linear  to  quadratic  to  cubic,  etc.  and  are 
called  hierarchical  because  their  contributions  to  the  accuracy  of 
the  solution  will  be  of  diminishing  importance  (Reference  7). 

(3)  use  both  h  and  p  refinements  simultaneously. 

The  h  and  p  methods  for  improving  finite  element  approximations  are 
compared  schematically  in  Figure  5. 

3.  The  Hierarchical  Approach 

The  hierarchical  approach  of  finite  element  analysis  has  been  used  pre¬ 
viously  to  seek  converged  solutions  for  structural  mechanics  problems  (Ref¬ 
erences  8-10).  Comparisons  of  finite  element  convergence  via  the  h 
approach  (mesh  refinement)  and  via  the  p  approach  (increasing  polynomial 
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ructural  solution 

sequence  of  analyses  (if  needed) 

accuracy  improved  via  structural  hierarchical  elements 

.  Hierarchical  integrated  thermal -structure  analysis. 


order)  are  made  in  Reference  9,  and  the  p  approach  is  found  to  be 
superior.  Reference  10  advocates  the  general  adoption  of  the  hierarchical 
approach  and  addresses  questions  of  adaptive  refinement  and  error  estima¬ 
tion. 

The  hierarchical  finite  element  approach  for  integrated  thermal  struc¬ 
tural  analysis  uses  a  common  discretization  for  the  analyses  and  seeks 
improvement  in  the  effectiveness  of  the  analyses  by:  (1)  improving  the 
accuracy  of  the  thermal  analysis  by  using  hierarchical  temperature  interpol¬ 
ation  functions  to  converge  the  thermal  solution,  (2)  using  the  converged 
temperature  distribution  to  compute  improved  finite  element  equivalent  nodal 
forces,  (3)  using  hierarchical  displacement  functions  (not  necessarily  the 
same  order  as  the  temperature  interpolation  functions)  to  converge  the 
structural  solution.  The  hierarchical  integrated  thermal-structural  analy¬ 
sis  is  shown  schematically  in  Figure  6. 

In  the  thermal  analysis,  element  temperatures  are  taken  in  the  form 

PT 

T(x,y,t)  =  Z  [  Nj(x,y)  ]  {T1  (t) }  (H) 

i=l 

where  [Nj(x,y}]  denote  temperature  interpolation  functions,  and  { T1 ( t ) J 
denote  vectors  of  unknown  temperatures.  For  i  =  1,  the  first  term  on  the 
right  hand  side  of  Equation  14  includes  the  first  order  bilinear  interpol¬ 
ation  functions,  and  the  vector  {T1}  denotes  the  nodal  temperatures.  For 
i  >  2  the  next  terms  represent  the  hierarchical  interpolation  where  the 
superscript  i  denotes  the  order  of  the  interpolating  polynomial  varying 
from  i  =  2  to  the  highest  order  for  the  thermal  analysis,  Pjt  and  the 
vectors  { T1 }  for  i  >  2  denote  the  nodeless  variables  associated  with  each 
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Figure  6.  Convergence  of  finite  element  approximations. 


interpolating  polynomial.  Typical  hierarchical  temperature  interpolating 
polynomials  and  nodeless  variables  are  shown  in  Figure  7. 

In  the  structural  analysis,  element  displacements  are  taken  in  the  form 

PS 

{u(x,y,t}}  =  I  CnZ  (x,y)  ]{  u1  ( t) }  (15) 

1=1  b 

where  [N^(x,y)]  denote  displacement  interpolation  functions,  and  { u 1 ( t ) } 
denote  vectors  of  unknown  displacements.  If  i  =  1,  the  first  term  on  the 
right  hand  side  of  Equation  15  denotes  the  first  order  interpolation  func¬ 
tions;  the  next  terms  represent  the  hierarchical  displacement  interpolation 
where  i  denotes  the  order  of  interpolating  polynomial.  For  i  >  2,  the 
vectors  {  u1 { t) }  represent  the  nodeless  structural  variables.  Note  that 
the  highest  order  P<.  for  the  structural  analysis  may  differ  from  the  high¬ 
est  order  in  the  thermal  analysis  Py  to  permit  optimization  of  each  solu¬ 
tion  independently.  Typical  hierarchical  displacement  interpolation  poly¬ 
nomials  and  nodeless  variables  are  shown  in  Figure  8. 

4.  Comments  on  the  Approach 

The  hierarchical  approach  for  integrating  thermal-structural  analysis 
is  quite  general  and  has  significant  potential  for  further  development.  The 
approach  is  based  on  using  well-known  finite  element  methodology  principally 
from  structural  analysis.  For  example,  the  interpolation  functions  shown  in 
Fiqures  7-8  are  the  well  known  Serendipity  functions  (Reference  6).  How¬ 
ever,  a  significant  effort  is  required  for  the  development  of  optimum  com¬ 
puter  programs  to  implement  the  concept.  In  this  study  only  very  primitive 
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VARIABLE  TYPICAL  INTERPOLATION  FUNCTIONS 
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Figure  7.  Typical  hierarchical  interpolation  functions  for  thermal  analysis. 
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Figure  8.  Typical  hierarchical  interpolation  functions  for  structural  analysis. 


computer  codes  were  developed  for  the  numerical  applications.  No  effort  was 
made  to  develop  efficient  general  codes  because  of  the  modest  funding  of 
this  study.  Instead  exploratory  programs  were  written  to  demonstrate  the 
concept.  For  instance,  the  highest  order  interpolation  functions  implement¬ 
ed  was  quadratic. 
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SECTION  IV 


NODE LESS  PARAMETER  STRUCTURAL  FINITE  ELEMENTS 
The  development  of  a  new  structural  element  using  the  nodeless  param¬ 
eter  approach  is  presented  in  this  section.  The  nodeless  parameter  approach 
can  reduce  thermal  stress  discontinuities  between  elements  that  are  produced 
by  conventional  finite  elements  as  demonstrated  in  Section  1 1.3.  The 
approach  improves  stress  distributions  without  adding  additional  element 
unknowns.  The  formulation  for  a  rectangular  element  with  linear  temperature 
distribution  is  first  presented.  The  approach  is  then  extended  to  quadri¬ 
lateral  and  hexahedral  elements.  Finally,  the  formulation  for  a  rectangular 
element  with  quadratic  temperature  distribution  (PT  =  2)  is  derived  to 
demonstrate  the  approach  for  hierarchical  thermal  element  temperature 
interpolation  polynomials. 

1.  Rectangular  Element  Formulation  for  Linear  Temperature  Distribution 
1.1  Element  Interpolation  Functions 

For  the  one-dimensional  nodeless  parameter  structural  element  presented 
in  Reference  5,  the  element  displacement  distribution  was  derived  from  the 
governing  differential  equation  of  a  thermo-elastic  member.  The  displace¬ 
ment  {$}  is  written  as  the  combination  of  the  conventional  linear  dis¬ 
placement  distribution  and  a  displacement  distribution  associated  with  ele¬ 
ment  nodal  temperatures, 

1*1  =  CNs1{6}  +  [N]{T  -  Tref}  (16) 

where  [N^]  and  (6)  are  the  conventional  element  displacement  interpo- 
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1  at i on  matrix  and  the  element  nodal  displacement  vector,  respectively.  In 
the  second  term,  [N]  represents  the  nodeless  parameter  interpolation 
matrix,  and  {t}  is  the  vector  of  element  nodal  temperatures.  The  nodal 
temperatures  are  known  and  have  the  role  of  parameters  in  the  structural 
analysis.  For  two-  and  three-dimensional  problems,  it  is  not  possible  to 
derive  the  element  displacement  distribution  in  a  simple  closed  form  solu¬ 
tion  from  the  equations  of  elasticity.  However,  the  idea  of  assuming  the 
element  displacement  distribution  in  the  form  of  Equation  16  can  be  used 
where  the  nodeless  parameter  interpolation  matrix  [N]  is  to  be  determined 
by  an  alternate  approach. 

For  a  two-dimensional  four-node  structural  element,  the  element  dis¬ 
placement  distributions  in  the  form  of  Equation  16  are. 
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where  u  and  v  are  the  elenent  in-plane  displacement  distributions  in 
local  coordinates  x  and  y. 

For  a  rectangular  element,  the  first  term  on  the  right  hand  side  of 
Equation  17  represents  the  bilinear  displacement  distributions  which  are 
given  by  Equation  6  and  are  illustrated  in  Figure  3.  With  the  assumed  dis¬ 
placement  distribution  in  the  form  of  Equation  16,  the  element  strain-dis¬ 
placement  relations  are 

{<=}  *  [Bs]  {6}  +  [B]  {T  -  Tref}  (18) 

where  [B]  is  the  nodeless  parameter  strai n-di spl acement  matrix  that  is  to 
be  determined.  Using  the  thermo-elastic  stress-strain  relations,  the  ele¬ 
ment  stresses  are, 

("1  -  [D]  { [B$]  {«}  +  [B]  {T  -  Tref}  -  M  (T(x.y)  -  Tref)J  (19) 

In  the  case  of  the  plane  stress,  for  example,  the  above  element  stresses  can 
be  written  in  terms  of  strains  as. 
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where  e  ,  e  and  y  are  the  strain  components  corresponding  to  the 
x  y  xy 

additional  nodeless  parameter  interpolations  introduced  in  Equation  16. 

As  mentioned  in  section  II. 3,  the  discontinuity  of  stress  a x  between 
elements  is  produced  because  the  thermal  strain  aT  is  one  order  higher 
algebraically  than  the  strain  ex .  The  stress  discontinuity  can  be  reduced 
by  forcing  the  additional  strain  term  ex  to  have  the  same  algebraic  order 
as  the  thermal  strain  aT,  i.e. 


(x,y)  -  a(T(x,y)  -  Tref)  +  C  (T(x,y)  -  Tr„f) 


ref  ‘ 


(21) 


where  C  is  a  constant.  This  expression  can  be  written  explicitly  as. 
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where  ,  i  =  1,4  are  the  nodeless  parameter  u-displ acement  interpola¬ 
tion  functions  shown  in  Equation  17,  These  interpolation  functions  are 
unknown;  however,  they  may  be  determined  by  applying  the  method  of  undeter¬ 
mined  coefficients  to  the  above  equation.  For  example,  by  comparing  the 
coefficients  of  TA : 


—  =  a  (1  -  JL)  (1  -  JL)  +  a, 

9x  a  b 

a  solution  of  the  nodeless  parameter  interpolation  function  is, 

2 

N  i  =  a  (x  -  — )  (1  — +  a*x  +  a2 
U1  2a  b 

where  ai  and  a2  are  constants  to  be  determined  by  requiring  that  the 
displacement  u  be  continuous  between  elements.  To  satisfy  this  continuity 
requirement,  the  above  nodeless  parameter  interpolation  function  must  vanish 
at  x  =  0  and  x  =  a,  i  ,e . 


Nul(x  =  0, y)  =  0 
%i(x  =  a,y)  =  0 


(23) 


With  these  two  boundary  conditions  the  constants  ai  and  a2  can  be  deter¬ 
mined  and  the  nodeless  parameter  interpolation  function  is  obtained  as. 
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Similarly,  the  unknown  nodeless  parameter  interpolation  functions  N^,  Nu3 
and  can  be  derived  using  the  boundary  conditions  given  in  Equation  23. 

Therefore,  the  element  u-di spl acement  distribution  shown  in  Equation  17 
becomes 
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For  an  isothermal  problem,  the  element  u-di spl acement  distribution 
shown  above  reduces  to  the  bilinear  distribution  as  assumed  in  the  conven¬ 
tional  element.  A  comparison  of  the  u-di spl acement  distributions  for  a 
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bilinear  conventional  element  (Equation  6)  and  a  nodeless  parameter  element 
(Equation  25)  is  shown  in  Figure  9(a). 

The  nodeless  parameter  interpolation  functions  for  the  element  v-dis- 
pl acement,  i  =  1,4,  can  be  derived  using  the  same  procedure.  The 

appropriate  boundary  conditions  which  preserve  the  element  v-displ  acement 
continuity  are. 


Nvi  (x,y=0)  =  0 
Nvi  (x,y=b)  =  0 


(26) 


The  element  v-displacement  distribution  with  nodeless  parameters  is 


v(x,y) 


^)(1  -  -*) 
a  b 


JL  (1  -  J)  AJL 
a  b  a  b 
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-  Conventional  bilinear  displacements 

-  Nodeless  parameter  displacements  for 
bilinear  temperatures 
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Figure  9.  Rectangular  structural  element  displacement  distributions. 


The  elenent  v-displ acement  distribution  is  compared  with  the  bilinear  con¬ 
ventional  element  distribution  in  Figure  9(b). 

1 . 2  Element  Stiffness  Matrix  and  Thermal  Force  Vectors 

Finite  element  equations  and  element  matrices  can  be  derived  by  apply¬ 
ing  the  principle  of  minimum  potential  energy  (Reference  6).  For  simplicity 
in  the  derivation  here,  the  effects  of  external  applied  forces  are  not  in¬ 
cluded;  only  thermal  effects  are  considered. 

The  total  potential  energy  is  represented  by  the  internal  strain  energy 
written  in  the  form. 


1  /  [t  -  «0]  {c}  dV 

2  V 


(28) 


n 


where  V  is  the  element  volume,  and  { a }  denotes  the  vector  of  stress  com¬ 
ponents;  [e]  denotes  the  row  matrix  of  the  total  strain  components  shown 

in  Equation  18,  and  [eo]  denotes  the  row  matrix  of  the  thermal  strain  com¬ 
ponents.  For  example,  the  vector  of  thermal  strain  components  for  the  plane 
stress  problem  is, 


*  <T  -  Tref> 

{e0}  »  1°!  (T  -  Tref)  -  (T  -  Tref) 


(29) 


0 


Using  the  stress-strain  relations, 


{ o}  =  [D]  je  -  eol 


(30) 
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the  total  potential  energy.  Equation  28,  becomes 


n  =  ~~  !  [e  -  e0]  [D]  {e  -  c0J  dV  (31) 

8y  substituting  the  strain-displacement  relations  for  the  nodeless  parameter 
element.  Equation  18,  the  total  potential  energy  can  be  written  in  the  form, 

"  ■  -^(«lT  /„  [8s]T[D][Bs]dV  {«]  -{«}T/  [Bs]T[D]  {t0)  dV 

+  _i  [e0l[D]  (E0)  dV  +  {S)T  /  [Bs]T[D][B]  { T-Tref}dV)  (32) 

The  element  equilibrium  equations  are  then  derived  by  performing  the  mini¬ 
mization  of  the  total  potential  energy  with  respect  to  element  modal 
unknowns.  The  equilibrium  equations  have  the  form, 

[Ks]  {6}  =  {  Ft}  -  {Ft}  (33) 

where  [K<.]  and  {  F^l  are  the  element  stiffness  matrix  and  the  equivalent 
nodal  thermal  load  vector,  respectively.  These  two  element  matrices  are 
identical  to  those  obtained  from  the  conventional  element  formulation  given 
in  Equation  9a  and  9b,  respectively .  The  vector  { Fy}  is  produced  by  the 
element  nodeless  parameter  interpolation  functions  and  is  given  by 


|Ft1  -  /  [8s]T[0][B]  {T-Trefl  dV 


(34) 


As  shown  by  the  element  equi 1 ibriin  equations,  the  nodeless  parameter 
approach  does  not  require  extra  element  unknowns  or  modify  the  stiffness 
matrix,  however,  an  additional  element  nodal  force  vector  { Fj }  is  intro¬ 
duced.  Due  to  this  additional  element  nodal  force  vector,  nodal  displace¬ 
ments  computed  are  different  in  general  from  the  conventional  finite  element 
sol ution. 

1.3  Element  Stresses 

Once  element  nodal  displacements  {6}  are  obtained,  element  stresses 
can  be  computed  from  the  thermo-elastic  stress-strain  relations, 

(a)  =  [D]  {[Bs]  {&}  +  [B]  <T-Tref)  -  {a}  (T-Tref)}  (35) 

2.  Quadrilateral  Element  Formulation 

A  two-dimensional,  four-node  elenent  with  a  general  quadrilateral  shape 
is  shown  in  Figure  10(a).  In  general,  the  finite  element  matrices  are 
integrals  over  the  elenent  volume,  and  the  integrations  are  performed  using 
Gauss  quadrature  in  the  £-n  coordinates  shown  in  Figure  10(b).  The  ele¬ 
ment  interpolation  functions  are  defined  in  terms  of  the  natural  coordinate 
variables  £,  and  n.  Details  of  the  coordinate  transformation  for  conven¬ 
tional  interpol ation  functions  are  found  in  Reference  6. 

To  obtain  the  nodeless  parameter  interpolation  functions  in  terms  of 
the  natural  coordinates  ?  and  n,  the  nodeless  parameter  interpolation 
functions  for  the  rectangular  element  given  in  Equations  25  and  27  are  used. 
The  relations  between  the  local  x-y  coordinates  (Figure  3)  and  the  £,-n 
coordinates  (Figure  9b)  for  the  rectangular  element  are 
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Figure  10.  Four-node  isoparametric  finite  element  in  global 
and  natural  coordinates. 


36 
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(3  6a) 


and 
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~  (i  +  o 
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-l  ■  —  (1  +  n) 
b  2 


(36b) 


With  these  relations,  the  element  displacement  distributions  become, 


uU.n)  =  [Ni  N2  N3  M*] 


v(5,n)  =  [N1  Ifc  N3  N,] 
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Vi 
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Vef 
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(37a) 


(37b) 


where  N.,  i  =  1,4  are  the  conventional  bilinear  interpolation  functions. 


Ni  =  —  (1-  £ )  (1-ti) 
4 

N3  *  _L  (1+  5 )  (l+n) 
4 


Ng  =  —  (1+  O(l-Tl) 
4 

N4  =  —  (1-  C )  (1  +r>) 


(38a) 


and  Nui  and  Nv>,  i  =  1,4,  are  the  nodeless  parameter  interpolation 
functions  corresponding  to  the  u-  and  v-di spl acements,  respectively. 
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N  i  =  —  a(l-€2 ) (1-n ) 

ul  16 


N,  =-—  a(l-C2 ) (1+n) 
16 


—  b(lH)(l-n2) 

16 

(38c) 

—  b(l-5 )  (1-n2  ) 

16 


Nvl  “  Y&  bC1_0Cl'n2)  Nv2  = 

^v3  =  ““  bCl^Xl-n2)  iv4=- 


N?  =  -—  a(l-S2)(l-n) 

^  16 

(38b) 

Nu4  =  —  a(l-?2){l+n) 

16 


Note  that  these  nodeless  parameter  interpolation  functions  depend  on  the  co¬ 
efficient  of  thermal  expansion,  the  dimensions  of  the  element  and  the 
natural  coordinates  K  and  n. 

For  an  element  with  general  quadri 1 ateral  shape,  it  is  not  possible  to 
derive  the  nodeless  parameter  interpolation  functions  as  for  the  rectangular 
element  given  by  Equations  38b-c,  However,  the  nodeless  parameter  interpo¬ 
lation  functions  obtained  for  rectangular  element  can  be  modified  and  used 
approximately  for  the  quadri 1 ateral  element  in  the  following  forms. 


Nul"  77  *12  Cl-c2) (1-n) 
16 

Nu3  =  *34  U“62 )  (l+n) 

16 


Nu2  = 

1  0 

—  *34  (l-ez) (l*n) 
16 


{39a} 
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N„1  =  *4  1  Cl-c )  (1-n2 ) 

16 


Nv2  =  77  *23  d+^Hl-1*2) 

lb 


(39b) 


N  3  =  -  —  ^23  (l+£)U-n2) 
VJ  16 


nv4  =  -  ~  *4i  (l-e)U-n2) 

1 6 


i  and  j.  With  the  above  nodeless  interpolation  functions,  the  element 
displacement  distributions  given  by  Equation  37  are  shown  in  Figure  11. 
3.  Hexahedral  Element  Formulation 

Figure  12  shows  an  eight-node  hexahedral  element  in  both  global  and 
natural  coordinates.  In  the  thermal  analysis,  the  conventional  element 
temperature  distribution  (Py  =  1)  is  expressed  in  terms  of  the  natural 
coordinates  £,  n  and  c  in  the  form  (Reference  6), 


8 


T(e,n,C,t)  =  [NyU.n.c)]  |T(t)}  = 


(40) 


where  N.,  i=l,8,  are  the  element  temperature  interpolation  functions 

given  by 


N.  =  J.  (l+  551)(l+nni)(l+«i) 
8 


(41) 


and  T-  are  the  element  nodal  temperatures  which  may  be  a  function  of  time 
t. 

The  congruent  structural  element  at  each  node  has  three  displacement 
unknowns,  u,  v,  and  w  in  the  directions  of  the  element  local  coordin- 
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Conventional,  bilinear.  4-node  element 


40 


Figure  11.  Quadrilateral  structural  element  displacement  distributions. 
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(a)  Global  coordinates 


Figure  12.  Eight-node  isoparametric  finite  element  in  global 
and  natural  coordinates. 
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ates  x,  y  and  z,  respectively.  The  nodeless  parameter  structural  ele¬ 
ment  can  be  constructed  by  first  writing  the  three  displacement  distribu¬ 
tions  in  the  following  forms, 

u  =  [Ns]  {u>  +  [Nu]  {T-Tref}  (42a) 

v  =  [Ns]  {v}  +  [Nv]  {T-Tref>  (42b) 

and  w  =  [Ns]  {w}  +  [N^]  {T-Tref>  (42c) 

where  [N<-]  is  the  row  matrix  of  the  conventional  displacement  interpola¬ 
tion  functions  given  by  Equation  41.  [N  ],  [N  ]  and  [N w3  are  the  node¬ 
less  interpolation  matrices  to  be  determined,  and  {T}  is  the  vector  of 
element  nodal  temperatures  shown  in  Equation  40. 

To  determine  the  nodeless  interpolation  functions,  the  procedures 
described  in  sections  IV. 1  and  IV. 2  are  employed.  Using  the  assumed  element 
displacenent  distributions  in  the  form  of  Equation  42,  the  element  strain- 
displacement  and  stress-strain  relations  are  formulated.  To  reduce  the  ele¬ 
ment  stress  discontinuity,  the  strains  associated  with  the  nodeless  param¬ 
eter  interpolation  functions  are  forced  to  have  the  same  algebraic  order  as 
the  thermal  strains.  These  conditions  and  the  element  boundary  conditions 
required  to  preserve  the  element  displacement  continuity  lead  to  the  deter¬ 
mination  of  the  unknown  nodeless  parameter  interpolation  functions.  The 
nodeless  parameter  interpolation  functions  corresponding  to  the  u,  v  and 


42 


w  displ  acements  are: 


Nul  =  *  N,2 


Nu5  =  "  Nu6 


Nu7  "  '  Nu8 


=  —  *12  (l-C2)(l-n)(l-0 
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=  -_2  *34  (l-?2)(l«)(l-5) 
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215  (1-C)  (1-T!)  (1-C2  ) 

*26  (1+e) (i-n)(l-^2) 

*37  (lH)(l+n)(l-?2} 

*48  d-e)a+Ti)(i-«4) 


(4  3c) 


where  *  .  j , 
and  j . 


i  ,j  =  1,8,  are  the  length  of  element  edges  between  nodes  i 


Element  matrices  can  be  derived  as  described  in  section  IV. 1.2  and  the 
element  equations  are  obtained  in  the  form  of  Equation  33.  Since  each  node 
has  three  displacement  components,  there  are  24  unknowns  for  an  element 
which  is  equal  to  the  number  of  unknowns  used  in  the  conventional  8-node 
element.  After  the  nodal  displacement  components  are  computed,  the  element 
displacement  distributions  and  the  element  stresses  are  obtained  using 
Equations  42  and  35,  respectively. 

4.  Formulation  for  Quadratic  Temperature  Distribution 

In  the  thermal  analysis,  a  two-dimensional  element  with  quadratic  temp¬ 
erature  distribution  can  be  constructed  using  the  hierarchical  approach.  As 
described  in  Section  III. 3,  the  hierarchical  element  temperature  distribu¬ 
tion  (Pj.  =  2)  is  expressed  in  the  form. 
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T  =  [N,  N2  N3  M 


(44) 
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1  +  [N5  n6  n7  n8]  ( 

t6 

t3 

T? 

T4 

Ts 

where  T,. ,  i  =  1,4  are  the  element  nodal  temperatures  and  T.,  i  =  5,8 
are  the  nodeless  variables.  For  a  rectangular  element,  the  element  interpo¬ 
lation  functions  N. ,  i  =  1,4  are  the  same  as  for  the  conventional  bi¬ 
linear  four-node  element  given  in  Equation  2.  The  interpolation  functions 
N . ,  i  =  5,3  for  the  nodeless  variables  are  given  by. 


n6  -  4  JL  J.  (1  -  J) 

a  b  b 

(45) 

N8  =  4(1-—)  1  (1  -  JL) 
a  b  b 


To  formulate  a  nodeless  parameter  structural  element  for  the  quadratic 
temperature  distribution,  the  procedures  described  in  Section  IV. 1.1  are 
employed.  The  structural  element  displacement  distributions  are  first 
written  as  the  combination  of  the  bilinear  distribution  and  the  distribution 
associated  with  the  element  temperature. 
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v  =  [Ni  N2  N3  N,  ] 
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where  and  N  ^  i  =  1,3  are  the  nodeless  parameter  interpolation 

functions  to  be  determined.  The  strain-displacement  and  the  stress-strain 
relations  shown  in  Equations  18  and  19  are  then  formulated.  To  reduce  the 
stress  discontinuity  between  elements,  the  strains  produced  by  the  nodeless 
parameter  displacement  distributions  are  forced  to  have  the  same  order  as 
the  thermal  strains.  As  an  example  for  plane  stress,  this  requirement  for 
the  u  displacement  yields 
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where  the  coefficients  in  the  row  matrix  of  the  second  term  on  the  right 
hand  side  of  the  equation  are  constants.  Using  the  method  of  undetermined 
coefficients  and  imposing  the  conditions  of  displacement  continuity  require¬ 
ment  between  elements,  the  unknown  nodeless  parameter  interpolation  func¬ 
tions  *N  ^ ,  i  =  1,8  are  determined.  Tbe  nodeless  parameter  interpolation 
functions  corresponding  to  the  u  and  v  displacements  are. 
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5.  Comments  on  Formulation 

For  the  two-dimensional  element  described  in  Section  IV. 1,  the  nodeless 
parameter  interpolation  functions  were  obtained  for  a  rectangular  element. 
For  a  general  quadrilateral  element,  the  nodeless  parameter  interpolation 
functions  do  not  appear  to  be  derived  easily.  The  nodeless  parameter  inter¬ 
polation  functions  for  rectangular  element  were  modified  and  used  approxi¬ 
mately  for  a  quadrilateral  element.  For  this  reason,  the  quadrilateral 
nodeless  parameter  element  (Equations  37-39}  will  provide  maximum  perform¬ 
ance  if  the  element  shape  is  close  to  rectangular. 
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Figure  9  shows  that  the  nodeless  parameter  interpolation  functions  are 
not  geometrically  isotropic;  the  functions  depend  on  the  element  orienta¬ 
tion.  This  means  that  care  must  be  taken  in  element  nodal  numbering  to  en¬ 
sure  compatibility  of  displacements  at  element  interfaces.  This  lack  of 
geometric  isotropy  restricts  the  generality  of  the  element  formulation.  The 
element  displacement  and  stress  distributions  compared  to  the  conventional 
element  solutions  may  not  be  improved  when  this  situation  occurs  and  further 
investigation  is  required  to  evaluate  the  element  performance  under  these 
conditions. 
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SECTION  V 


APPLICATIONS 

To  demonstrate  the  capabilities  of  the  hierarchical  and  the  nodeless 
parameter  finite  elanents  described  in  Sections  III  and  IV,  four  thermal- 
stress  problems  are  analyzed:  (1)  a  free  expansion  plate  with  linear 
temperature  distribution,  (2)  a  fixed  end  beam  with  quadratic  temperature 
distribution  through  the  beam  depth,  (3}  a  simplified  wing  section  with 
aerodynamic  heating,  and  (4)  a  convectively  cooled  laser  mirror.  In  each 
problem,  benefits  of  utilizing  the  hierarchical  and  nodeless  parameter 
finite  elements  are  demonstrated  by  comparison  with  results  from 
conventional  finite  elements  and,  where  possible,  analytical  solutions. 

1.  Free  Expansion  Plate  with  Linear  Temperature  Distribution 

8oth  hierarchical  and  nodeless  parameter  finite  elements  are  used  for 
the  structural  analysis  of  the  free  expansion  plate  with  linear  temperature 
distribution  (Figure  4(a))  as  described  previously  in  Section  II. 3.  The 
plate  is  modeled  by:  (l)  one  conventional  finite  element,  (2)  one  nodeless 
parameter  finite  elanent,  and  (3)  one  hierarchical  finite  element  (P^  =  2). 
To  illustrate  the  performance  of  these  finite  elements,  the  analytical  solu¬ 
tions  for  displacements  and  stresses  were  derived.  Since  the  plate  is  sub¬ 
jected  to  a  linear  temperature  distribution  and  is  free  to  expand,  all 
stress  components  are  zero.  With  the  linear  temperature  distribution  in  the 
x-direction,  the  analytical  solution  for  the  displacement  are: 

u  =  a  x  TA  +  ~  (x2  -  y2)  (Tz-  Ti)  (48a) 
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and 


v  =  oy  Ti  +  —  xy  (T2  *  T 1 ) 
a 


(48b) 


where  a  is  the  thermal  expansion  coefficient  of  the  plate;  T*  and  T2 
are  the  temperatures  along  the  edges  x  s  0  and  x  =  a,  respectively. 

Figure  13  shows  the  comparative  u-displ acement  distributions  for  the 
analytical  solution  and  the  three  finite  element  solutions  along  the  edges 
y  =  0  and  y  =  b.  Since  the  hierarchical  element  with  P$*  2  uses 
quadratic  displacement  distributions,  the  displacement  solutions  obtained 
are  exact.  At  y  =  0  the  computed  nodal  u-di spl acenent  at  x  =  a  from  all 
three  finite  element  models  are  exact;  however,  the  conventional  finite 
element  is  unable  to  provide  a  realistic  displacement  distribution  for  0  < 
x  <  a.  Similarly,  along  edge  y  =  b,  the  conventional  finite  element 
solution  tends  to  average  the  true  displacement  distribution  whereas  the 
nodeless  parameter  element  gives  excellent  agreement  with  the  exact 
solution. 

El  ament  stress  distributions  a  along  the  edge  y  =  0  and  o  along 

x  y 

the  edge  x  =  0  obtained  from  the  three  finite  element  models  are  shown  in 
Figure  14(a)  and  (b),  respectively.  The  hierarchical  element  provides  exact 
zero  stresses  throughout  the  plate.  Due  to  the  linear  temperature  distri¬ 
bution  in  the  x-direction,  an  unrealistic  distribution  for  a  is  obtained 
from  the  conventional  element.  The  nodeless  parameter  element  improves  this 
stress  distribution  significantly  and  clearly  demonstrates  the  capability  of 
providing  more  realistic  stress  distributions. 

2.  Fixed  End  Beam  with  Quadratic  Temperature  Distribution 

A  beam  with  a  symmetrical  quadratic  temperature  distribution  through 
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temperature  distribution. 


V0-*1  -  ANALYTICAL  AND  HIERARCHICAL  F  E.  (P  =2) 
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Figure  14.  Comparative  thermal  stress  distributions  for  a  free  expansion  plate  with  linear 
temperature  distribution. 


its  depth  is  shown  in  Figure  15(a).  Hie  beam  is  fixed  in  the  longitudinal 
direction  at  both  ends  but  is  free  to  move  in  the  transverse  direction  due 
to  the  constraints  u  =  0  and  v  ~  v(y).  Using  synmetry,  the  upper  half  of 
the  beam  is  modeled  by:  (1)  one  conventional  element,  and  (2)  one  nodeless 
parameter  element. 

Fiqure  15(b)  shows  the  comparative  transverse  displacement  distribu¬ 
tions  from  an  analytical  solution  (Reference  11),  the  conventional  element 
and  the  nodeless  parameter  element  solutions.  The  conventional  element 
gives  a  relatively  high  error  for  both  the  nodal  displacement  and  the  dis¬ 
placement  distribution.  The  nodeless  parameter  element  yields  the  exact 
nodal  displacement  and  gives  an  excellent  representation  of  the  displacement 
distribution. 

The  stress  distributions  obtained  from  these  finite  elenent  models  are 
compared  with  the  analytical  solution  in  Figure  15(c).  The  element  stress 
distribution  predicted  by  the  nodeless  parameter  element  is  in  very  good 
agreement  with  the  analytical  solution.  The  conventional  element  is  unable 
to  provide  details  of  the  nonuniform  stress  distribution.  As  mentioned 
earlier,  element  centroidal  stresses  are  not  the  correct  estimates  of  the 
true  stresses  in  general.  This  phenomenon  is  clearly  shown  by  this  example 
where  the  centroidal  stress  obtained  from  the  conventional  elements  is 
underestimated  significantly. 

3.  Simplified  Wing  Section  with  Aerodynamic  Heating 

A  simplified  wing  section  (Figure  16a)  consisting  of  metallic  top  and 
bottom  skins  connected  by  corrugated  spars  is  subjected  to  synmetrical,  non- 
uniform  step-function  representations  of  time-dependent  aerodynamic  heating. 
The  finite  element  discretization  is  based  on  a  unit  length  of  the  wing  sec- 
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(a)  Fixed  end  beam  with  nonlinear  temperature  distribution. 


lb)  Comparative  transverse  displacement  distributions. 


Figure  15.  Conventional  and  nodeless  parameter  finite  element  solutions 
for  a  fixed  end  beam  with  nonlinear  temperature  distribution. 
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(c)  Comparative  thermal  stress  distributions. 
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Figure  15.  Concluded. 
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(a)  SIMPLIFIED  WING  SECTION  WITH  AERODYNAMIC 
HEATING 


t i on  and  models  a  one-dimensional  transverse  variation  of  temperature  and 
thermal  stress.  Three  finite  element  models  were  used  for  the  analyses. 

The  first  model  consists  of  seven  linear  elenents  (Py  »  P^  =  1);  two  ele¬ 
ments  each  for  the  top  and  bottom  skins  and  one  elenent  for  each  spar.  The 
second  model  used  a  refined  mesh  with  ten  linear  elenents  for  each  skin. 

The  third  model  is  identical  to  the  first  model  except  quadratic  hierarchi¬ 
cal  interpolation  functions  were  employed  for  the  thermal  analysis  (Py  =  2) 
and  linear  interpol ation  functions  (P^  *  1)  were  employed  for  the  struc¬ 
tural  analysis.  The  comparative  corresponding  skin  temperatures  at  t  -  150 
s.  are  shown  in  Figure  16b,  and  the  corresponding  skin  thermal  stresses  are 
shown  in  Figure  16c. 

For  both  temperature  and  stress,  the  hierarchical  approach  with  Py  =  2 
and  P<-  =  1  predicts  realistic  solutions  and  gives  good  agreement  with  the 
results  from  the  refined  mesh  of  linear  elements.  The  first  model  with 
linear  elements  is  unable  to  represent  the  non  uni  form  temperature  and  stress 
distributions.  But  with  the  same  discretization,  the  addition  of  the  quad¬ 
ratic  interpolation  functions  for  the  temperature  is  clearly  sufficient  to 
resolve  the  pertinent  details  of  both  the  thermal  and  structural  solutions. 
4,  Convectively  Cooled  Laser  Mirror 

Figure  17(a)  illustrates  a  design  concept  of  a  convectively  cooled 
laser  mirror.  The  mirror  is  subjected  to  high  uniform  specified  surface 
heating.  To  reduce  stresses  and  deformations,  the  mirror  is  convectively 
cooled  by  flow  through  multiple  coolant  passages.  Using  symmetry,  the 
finite  element  discretization  shown  in  Figure  17(b)  was  established  for  a 
preliminary  thermal-stress  analysis. 

The  thermal  analysis  was  simplified  by  using  a  two  dimensional  model  of 
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(b)  COMPARATIVE  TEMPERATURE  DISTRIBUTIONS  ALONG 
CHORDWISE  DIRECTION.  t=  150  sec. 
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(c)  COMPARATIVE  THERMAL  STRESS  DISTRIBUTIONS  ALONG 
CHORDWISE  DIRECTION,  t=150  sec. 
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(a)  CONVECTIVELY  COOLED  LASER  MIRROR 
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(b)  FINITE  ELEMENT  MODEL 


Figure  17.  Thermal -structural  analysis  of  convectively  cooled  laser  mirror 


a  mirror  cross-section  with  specified  coolant  temperature.  Steady-state 
conduction  with  convection  heat  transfer  to  the  coolant  was  analyzed.  The 
bottom  mirror  surface  was  assumed  to  be  adiabatic,  and  the  left  and  right 
edges  of  the  mirror  have  zero  heat  transfer  by  symmetry.  Conventional  bi¬ 
linear  thermal  elements  (P|  =  1)  were  used  in  the  analysis.  Temperature 
distributions  at  four  sections  through  the  mirror  are  shown  in  Figure  18. 
Since  the  temperatures  obtained  are  relatively  smooth,  the  conventional  bi¬ 
linear  elanent  model  for  the  thermal  analysis  was  sufficient  and  no  addi¬ 
tional  refinement  was  required. 

For  the  structural  analysis,  a  state  of  plane  strain  was  assumed.  Two 
finite  elanent  analyses  with  the  same  discretization  as  the  thermal  model 
were  performed.  Conventional  bilinear  elements  (P^  *  1)  and  hierarchical 
biquadratic  elements  (P^  =  2)  were  used  for  the  analyses.  Figure  19  shows 
comparative  displacement  distributions  along  the  mirror  surface.  The  hier¬ 
archical  elements  (Pc  =  2)  provide  a  smooth  displacement  distribution 
whereas  relatively  high  discontinuities  of  displacement  gradients  are  pro¬ 
duced  by  the  conventional  elements  (P<-  -  1).  Such  discontinuities  indicate 
the  need  for  mesh  refinement  if  conventional  elements  are  to  be  used. 

Thermal  stress  distributions  along  the  middle  plane  of  the  mirror  are 
shown  in  Figure  20.  The  conventional  elements  produce  element  stress  dis¬ 
continuities  whereas  a  more  realistic  stress  distribution  is  obtained  from 
the  hierarchical  elements.  The  thermal  stresses  from  the  hierarchical  ana¬ 
lysis  show  small  discontinuities  indicating  that  solution  convergance  is  not 
fully  attained  for  =  2,  suggesting  the  need  for  including  additional 
hierarchical  terms.  However,  the  principal  advantage  of  using  the  hier¬ 
archical  elements  is  clearly  demonstrated  by  this  example.  Enhanced  dis- 


63 


i 


64 


DISTRIBUTION  (b)  TRANSVERSE  TEMPERATURE  DISTRIBUTION 
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Figure  19.  Comparative  transverse  displacement  distributions 
along  mirror  surface. 
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Figure  20.  Comparative  mirror  in-plane  stress  distributions. 
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placement  and  thermal-stress  solutions  were  obtained  for  a  common  geometric 
model . 
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SECTION  VI 


CONCLUDING  REMARKS 

An  integrated  approach  for  enhanced  thermal-structural  analysis  was 
presented.  The  approach  focuses  on  applications  where  thermal  and  structur¬ 
al  models  often  differ  because  of  different  analysis  requirements.  The 
objectives  of  the  approach  are  to  provide  more  efficient  coupling  between 
the  thermal  and  structural  analysis  and  to  improve  the  accuracy  of  each  ana¬ 
lysis  particularly  the  thermal-stress  analysis.  The  integrated  approach  is 
based  on  usinq  the  same  geometric  model  with  a  common  nodal  discretization 
for  both  analysis  although  the  thermal  and  structural  models  can  employ  dif¬ 
ferent  elenents  to  suit  their  different  requirements. 

Two  approaches  for  integrating  finite  element  thermal  and  structural 
analyses  are  presented.  The  first  approach  is  based  on  applying  the  hier¬ 
archical  concept  of  finite  element  approximation  to  both  the  thermal  and 
structural  analysis.  In  an  hierarchical  approach  the  accuracy  of  the  finite 
element  approximation  is  improved  for  the  same  mesh  by  increasing  the  order 
of  interpolating  functions  and  introducing  additional  unknowns  via  nodeless 
variables.  The  hierarchical  approach  to  integrated  thermal-structural  ana¬ 
lysis  uses  a  common  discretization  for  the  thermal  and  structural  analysis 
and  seeks  improvements  in  the  accuracy  of  the  thermal  and  structural  analy¬ 
ses  by  independently  refining  the  solutions  using  hierarchical  interpol ation 
functions  for  successive  analyses.  A  key  step  in  coupling  the  analyses  is 
to  use  the  converged  temperature  distribution  to  compute  the  finite  elanent 
equivalent  thermal  forces.  The  second  approach,  called  a  nodeless  param¬ 
eter  approach,  uses  a  common  discretization  for  both  analyses  and  uses  hier¬ 
archical  interpolation  functions  to  converge  the  thermal  solution.  The 
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structural  analysis  is  based  on  using  new  temperature-dependent  displacement 
interpolation  functions  that  have  element  temperatures  as  parameters. 
Interpolation  functions  are  presented  for  a  rectangular  element,  and  the 
generalization  for  a  quadrilateral  is  briefly  discussed.  The  use  of  the  new 
interpolation  functions  can  improve  the  accuracy  of  the  structural  analysis 
without  adding  extra  unknowns. 

Four  two-dimensional  examples  are  presented  to  illustrate  the  two 
approaches.  Rectangular  elements  with  bilinear  interpolation  functions  are 
used  for  the  initial  analyses,  and  refinements  of  the  analyses  are  made 
using  the  hierarchical  approach  or  the  nodeless  parameter  approach.  Im¬ 
provements  in  the  accuracy  of  temperature  and  thermal-stresses  were  demon¬ 
strated  in  all  examples. 

The  hierarchical  approach  offers  the  greatest  potential  for  developing 
a  general  integrated  thermal -structural  analysis  method.  Maximum  flexi¬ 
bility  is  permitted  for  independently  improving  the  finite  element  approxi¬ 
mation  for  each  analysis  while  maintaining  a  common  discretization  and  the 
analyses  can  be  consistently  coupled  through  the  equivalent  thermal  forces. 
Additional  study  is  needed  to:  (1)  gain  experience  with  higher  order  inter¬ 
polation  functions,  (2)  develop  error  estimation  techniques  to  quantify 
conveyance,  and  (3)  study  computer  implementation  techniques.  The  nodeless 
parameter  approach  offers  the  advantage  of  improving  the  accuracy  of  the 
structural  analysis  without  adding  unknowns.  The  approach,  however,  needs 
additional  development  before  it  can  be  implemented  for  general  elements. 

The  integrated  thermal-stress  analysis  approach  based  on  hierarchical 
elements  provides  capability  to  improve  accuracy  and  efficiency  of  thermal- 
stress  analysis  of  complex  structures.  The  examples  presented  in  this 
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report  validate  some  basic  features  of  the  approach  for  two-dimensional 
thermal  stress  problems,  but  additional  research  is  needed  to  develop  the 
approach  to  its  full  potential. 
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LIST  OF  SYMBOLS 


a 

Surface  absorptivity 

[B] 

Nodeless  parameter  strain-displacement  interpolation 
matrix 

[Bs] 

Strain-displacement  interpolation  matrix 

[By] 

Temperature  gradient  interpolation  matrix 

c 

Specific  heat 

[C] 

Finite  element  capacitance  matrix 

[D] 

Elasticity  matrix 

E 

Modulus  of  elasticity 

ift} 

Finite  element  nodal  thermal  force  vector 

ifti 

Finite  element  nodeless  parameter  thermal  force 
vector 

[k] 

Thermal  conductivity  matrix 

[kc] 

Finite  element  conduction  matrix 

[Kr] 

Finite  element  radiation  matrix 

[KS] 

Finite  element  stiffness  matrix 

[N] 

Finite  element  interpolation  function  matrix 

m 

Finite  element  nodeless  parameter  interpolation 
function  matrix 

["Si 

Finite  element  displacement  interpolation  function 
matrix 

[Nr] 

Finite  element  temperature  interpolation  function 
matrix 

PT 

Order  of  interpolating  polynomial  for  thermal 
analysi  s 

PS 

Order  of  interpolating  polynomial  for  structural 
analysi  s 

q 

Surface  heating  rate 

7? 


*V 

Surface  incident  radiation  heating  rate 

{V 

Finite  element  surface  heating  load  vector 

{ Qr) 

Finite  element  radiation  heating  load  vector 

t 

Time 

*  T 

Temperature 

*  ^"ref 

Reference  temperature  for  zero  stress 

U,  V,  w 

Displacement  components 

V 

Vol  ume 

X,  y,  z 

Cartesian  coordinates 

{a} 

Vector  of  thermal  expansion  coefficients 

{6} 

Vector  of  finite  element  nodal  displacements 

P 

Mass  density 

a 

Stefan-Boltzmann  constant 

e 

Surface  emissivity 

it) 

Vector  of  finite  element  total  strains 

{eo^ 

Vector  of  finite  element  thermal  strains 

(a) 

Vector  of  finite  element  stresses 

V 

Poisson's  ratio 

n 

Internal  strain  energy 

£  ,  n  ,  c 

Finite  element  natural  coordinates 

s 
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